GTEx OTD joint heritability. Local h2 is estimated with SNPs within 1 Mb of each gene. distal h2 is estimated with SNPs that are eQTLs in the Framingham Heart Study on other chromosomes (FDR < 0.05).

  library(ggplot2)
  library(reshape2)
  library(dplyr)
  library(tidyr)
  library(GGally)
  library(grid)
  "%&%" = function(a,b) paste(a,b,sep="")
  source('/Volumes/im-lab/nas40t2/hwheeler/PrediXcan-Paper/scripts/Heather/make-figures/multiplot.R')
  my.dir <- '/Volumes/im-lab/nas40t2/hwheeler/cross-tissue/'
  fig.dir <- '~/GitHub/GenArch/GenArchPaper/Figures/'

Tissue-Specific joint local estimates

tislist <- scan(my.dir %&% 'TS.ten.tissue.list',sep="\n",what="character")
ts <- data.frame()
for(tis in tislist){
  data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")  
  loc.jt <- select(data,tissue,loc.jt.h2,loc.jt.se) %>% arrange(loc.jt.h2) %>% mutate(loc.jt.h2=ifelse(is.na(loc.jt.h2),0,loc.jt.h2)) %>% mutate(ymin = pmax(0, loc.jt.h2 - 2 * loc.jt.se), ymax = pmin(1, loc.jt.h2 + 2 * loc.jt.se)) %>% mutate(nonzeroCI=ymin>0,place=1:length(data$tissue))
  ts <- rbind(ts,loc.jt)
}

tsrmNA <- ts[complete.cases(ts),]

p<-ggplot(tsrmNA,aes(x=place,y=loc.jt.h2,ymin=ymin,ymax=ymax,color=nonzeroCI) ) + facet_wrap(~tissue,ncol=2) + geom_pointrange(col='gray')+geom_point()+ ylab(expression("local h"^2)) + xlab(expression("genes ordered by joint local h"^2))+theme_bw() + coord_cartesian(ylim=c(0,1)) + ggtitle("CT and TS")

###calc % nonzero for each tissue
### ADD to plot legend
pvec<-vector()
h2vec<-vector()
a<-tsrmNA %>% select(tissue,nonzeroCI) %>% spread(tissue,nonzeroCI)
for(i in 1:10){
  tis<-colnames(a)[i]
  cat("\n\n---",tis,"---\n")
  print(table(a[,i]))
  per <- signif(table(a[,i])/sum(table(a[,i])),3)*100
  print(per)
  pvec <- c(pvec,per[2])
}
## 
## 
## --- cross-tissue ---
## 
## FALSE  TRUE 
## 10442  2014 
## 
## FALSE  TRUE 
##  83.8  16.2 
## 
## 
## --- Adipose - Subcutaneous ---
## 
## FALSE  TRUE 
##  8983   100 
## 
## FALSE  TRUE 
##  98.9   1.1 
## 
## 
## --- Artery - Tibial ---
## 
## FALSE  TRUE 
## 10021   143 
## 
## FALSE  TRUE 
## 98.60  1.41 
## 
## 
## --- Heart - Left Ventricle ---
## 
## FALSE  TRUE 
##  9467   124 
## 
## FALSE  TRUE 
## 98.70  1.29 
## 
## 
## --- Lung ---
## 
## FALSE  TRUE 
##  9372    98 
## 
## FALSE  TRUE 
## 99.00  1.03 
## 
## 
## --- Muscle - Skeletal ---
## 
## FALSE  TRUE 
##  9820   211 
## 
## FALSE  TRUE 
##  97.9   2.1 
## 
## 
## --- Nerve - Tibial ---
## 
## FALSE  TRUE 
##  9665   184 
## 
## FALSE  TRUE 
## 98.10  1.87 
## 
## 
## --- Skin - Sun Exposed (Lower leg) ---
## 
## FALSE  TRUE 
## 11251   219 
## 
## FALSE  TRUE 
## 98.10  1.91 
## 
## 
## --- Thyroid ---
## 
## FALSE  TRUE 
##  9637   206 
## 
## FALSE  TRUE 
## 97.90  2.09 
## 
## 
## --- Whole Blood ---
## 
## FALSE  TRUE 
## 10764   244 
## 
## FALSE  TRUE 
## 97.80  2.22
###calc mean h2 for each tissue
a<-tsrmNA %>% select(tissue,loc.jt.h2) %>% spread(tissue,loc.jt.h2)
for(i in 1:10){
  tis<-colnames(a)[i]
  meanh2 <- signif(mean(a[,i],na.rm=TRUE),3)
  cat(tis,"mean h2:",meanh2,"\n")
  h2vec <- c(h2vec,meanh2)
}
## cross-tissue mean h2: 0.0553 
## Adipose - Subcutaneous mean h2: 0.0238 
## Artery - Tibial mean h2: 0.0251 
## Heart - Left Ventricle mean h2: 0.0319 
## Lung mean h2: 0.024 
## Muscle - Skeletal mean h2: 0.0229 
## Nerve - Tibial mean h2: 0.0292 
## Skin - Sun Exposed (Lower leg) mean h2: 0.0257 
## Thyroid mean h2: 0.0273 
## Whole Blood mean h2: 0.0237
ann_text <- data.frame( loc.jt.h2 = rep(0.78,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2<-p+geom_text(data=ann_text,aes(label=paste("mean_h^2 ==",mean_h2,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3) 
ann_text <- data.frame( loc.jt.h2 = rep(0.9,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2+geom_text(data=ann_text,aes(label=paste("percent_TRUE ==",percent,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)+ theme(legend.justification=c(0,1), legend.position=c(0,1))

Tissue-Specific joint global estimates

tislist <- scan(my.dir %&% 'TS.ten.tissue.list',sep="\n",what="character")
ts <- data.frame()
for(tis in tislist){
  data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")  
  glo.jt <- select(data,tissue,glo.jt.h2,glo.jt.se) %>% arrange(glo.jt.h2) %>% mutate(glo.jt.h2=ifelse(is.na(glo.jt.h2),0,glo.jt.h2)) %>% mutate(ymin = pmax(0, glo.jt.h2 - 2 * glo.jt.se), ymax = pmin(1, glo.jt.h2 + 2 * glo.jt.se)) %>% mutate(nonzeroCI=ymin>0,place=1:length(data$tissue))
  #glo.jt <- glo.jt[complete.cases(glo.jt),]
  ts <- rbind(ts,glo.jt)
}

tsrmNA <- ts[complete.cases(ts),]

p<-ggplot(tsrmNA,aes(x=place,y=glo.jt.h2,ymin=ymin,ymax=ymax,color=nonzeroCI) ) + facet_wrap(~tissue,ncol=2) + geom_pointrange(col='gray')+geom_point()+ylab(expression("global h"^2)) + xlab(expression("genes ordered by joint global h"^2))+theme_bw() + coord_cartesian(ylim=c(0,1)) + ggtitle("CT and TS")

###calc % nonzero for each tissue
### ADD to plot legend
pvec<-vector()
h2vec<-vector()
a<-tsrmNA %>% select(tissue,nonzeroCI) %>% spread(tissue,nonzeroCI)
for(i in 1:10){
  tis<-colnames(a)[i]
  cat("\n\n---",tis,"---\n")
  print(table(a[,i]))
  per <- signif(table(a[,i])/sum(table(a[,i])),3)*100
  print(per)
  pvec <- c(pvec,per[2])
}
## 
## 
## --- cross-tissue ---
## 
## FALSE  TRUE 
## 12123   333 
## 
## FALSE  TRUE 
## 97.30  2.67 
## 
## 
## --- Adipose - Subcutaneous ---
## 
## FALSE  TRUE 
##  9046    37 
## 
##  FALSE   TRUE 
## 99.600  0.407 
## 
## 
## --- Artery - Tibial ---
## 
## FALSE  TRUE 
## 10150    14 
## 
##  FALSE   TRUE 
## 99.900  0.138 
## 
## 
## --- Heart - Left Ventricle ---
## 
## FALSE 
##  9591 
## 
## FALSE 
##   100 
## 
## 
## --- Lung ---
## 
## FALSE  TRUE 
##  9463     7 
## 
##   FALSE    TRUE 
## 99.9000  0.0739 
## 
## 
## --- Muscle - Skeletal ---
## 
## FALSE  TRUE 
##  9844   187 
## 
## FALSE  TRUE 
## 98.10  1.86 
## 
## 
## --- Nerve - Tibial ---
## 
## FALSE  TRUE 
##  9848     1 
## 
##    FALSE     TRUE 
## 100.0000   0.0102 
## 
## 
## --- Skin - Sun Exposed (Lower leg) ---
## 
## FALSE  TRUE 
## 11355   115 
## 
## FALSE  TRUE 
##    99     1 
## 
## 
## --- Thyroid ---
## 
## FALSE  TRUE 
##  9831    12 
## 
##  FALSE   TRUE 
## 99.900  0.122 
## 
## 
## --- Whole Blood ---
## 
## FALSE  TRUE 
## 10882   126 
## 
## FALSE  TRUE 
## 98.90  1.14
###calc mean h2 for each tissue
a<-tsrmNA %>% select(tissue,glo.jt.h2) %>% spread(tissue,glo.jt.h2)
for(i in 1:10){
  tis<-colnames(a)[i]
  meanh2 <- signif(mean(a[,i],na.rm=TRUE),3)
  cat(tis,"mean h2:",meanh2,"\n")
  h2vec <- c(h2vec,meanh2)
}
## cross-tissue mean h2: 0.145 
## Adipose - Subcutaneous mean h2: 0.172 
## Artery - Tibial mean h2: 0.196 
## Heart - Left Ventricle mean h2: 0.287 
## Lung mean h2: 0.192 
## Muscle - Skeletal mean h2: 0.151 
## Nerve - Tibial mean h2: 0.257 
## Skin - Sun Exposed (Lower leg) mean h2: 0.224 
## Thyroid mean h2: 0.19 
## Whole Blood mean h2: 0.166
ann_text <- data.frame( glo.jt.h2 = rep(0.78,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2<-p+geom_text(data=ann_text,aes(label=paste("mean_h^2 ==",mean_h2,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3) 
ann_text <- data.frame( glo.jt.h2 = rep(0.9,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2+geom_text(data=ann_text,aes(label=paste("percent_TRUE ==",percent,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)+ theme(legend.justification=c(0,1), legend.position=c(0,1))

Tissue-Wide joint local estimates

tislist <- scan(my.dir %&% 'TW.ten.tissue.list',sep="\n",what="character")
tw <- data.frame()
for(tis in tislist){
  data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")  
  loc.jt <- select(data,tissue,loc.jt.h2,loc.jt.se) %>% arrange(loc.jt.h2) %>% mutate(loc.jt.h2=ifelse(is.na(loc.jt.h2),0,loc.jt.h2)) %>% mutate(ymin = pmax(0, loc.jt.h2 - 2 * loc.jt.se), ymax = pmin(1, loc.jt.h2 + 2 * loc.jt.se)) %>% mutate(nonzeroCI=ymin>0,place=1:length(data$tissue))
  #loc.jt <- loc.jt[complete.cases(loc.jt),]
  tw <- rbind(tw,loc.jt)
}

twrmNA <- tw[complete.cases(tw),]

p<-ggplot(twrmNA,aes(x=place,y=loc.jt.h2,ymin=ymin,ymax=ymax,color=nonzeroCI) ) + facet_wrap(~tissue,ncol=2) + geom_pointrange(col='gray')+geom_point()+ylab(expression("local h"^2)) + xlab(expression("genes ordered by joint local h"^2))+theme_bw() + coord_cartesian(ylim=c(0,1)) + ggtitle("CT and TW")

###calc % nonzero for each tissue
### ADD to plot legend
pvec<-vector()
h2vec<-vector()
a<-twrmNA %>% select(tissue,nonzeroCI) %>% spread(tissue,nonzeroCI)
for(i in 1:10){
  tis<-colnames(a)[i]
  cat("\n\n---",tis,"---\n")
  print(table(a[,i]))
  per <- signif(table(a[,i])/sum(table(a[,i])),3)*100
  print(per)
  pvec <- c(pvec,per[2])
}
## 
## 
## --- cross-tissue ---
## 
## FALSE  TRUE 
## 10442  2014 
## 
## FALSE  TRUE 
##  83.8  16.2 
## 
## 
## --- Adipose-Subcutaneous ---
## 
## FALSE  TRUE 
##  2981   588 
## 
## FALSE  TRUE 
##  83.5  16.5 
## 
## 
## --- Artery-Tibial ---
## 
## FALSE  TRUE 
##  2829   674 
## 
## FALSE  TRUE 
##  80.8  19.2 
## 
## 
## --- Heart-LeftVentricle ---
## 
## FALSE  TRUE 
##  1756   411 
## 
## FALSE  TRUE 
##    81    19 
## 
## 
## --- Lung ---
## 
## FALSE  TRUE 
##  3706   523 
## 
## FALSE  TRUE 
##  87.6  12.4 
## 
## 
## --- Muscle-Skeletal ---
## 
## FALSE  TRUE 
##  3450   606 
## 
## FALSE  TRUE 
##  85.1  14.9 
## 
## 
## --- Nerve-Tibial ---
## 
## FALSE  TRUE 
##  2267   800 
## 
## FALSE  TRUE 
##  73.9  26.1 
## 
## 
## --- Skin-SunExposed(Lowerleg) ---
## 
## FALSE  TRUE 
##  2917   622 
## 
## FALSE  TRUE 
##  82.4  17.6 
## 
## 
## --- Thyroid ---
## 
## FALSE  TRUE 
##  2786   724 
## 
## FALSE  TRUE 
##  79.4  20.6 
## 
## 
## --- WholeBlood ---
## 
## FALSE  TRUE 
##  3082   589 
## 
## FALSE  TRUE 
##    84    16
###calc mean h2 for each tissue
a<-twrmNA %>% select(tissue,loc.jt.h2) %>% spread(tissue,loc.jt.h2)
for(i in 1:10){
  tis<-colnames(a)[i]
  meanh2 <- signif(mean(a[,i],na.rm=TRUE),3)
  cat(tis,"mean h2:",meanh2,"\n")
  h2vec <- c(h2vec,meanh2)
}
## cross-tissue mean h2: 0.0553 
## Adipose-Subcutaneous mean h2: 0.0889 
## Artery-Tibial mean h2: 0.0978 
## Heart-LeftVentricle mean h2: 0.11 
## Lung mean h2: 0.0773 
## Muscle-Skeletal mean h2: 0.0715 
## Nerve-Tibial mean h2: 0.123 
## Skin-SunExposed(Lowerleg) mean h2: 0.0899 
## Thyroid mean h2: 0.105 
## WholeBlood mean h2: 0.0755
ann_text <- data.frame( loc.jt.h2 = rep(0.78,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2<-p+geom_text(data=ann_text,aes(label=paste("mean_h^2 ==",mean_h2,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3) 
ann_text <- data.frame( loc.jt.h2 = rep(0.9,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2+geom_text(data=ann_text,aes(label=paste("percent_TRUE ==",percent,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)+ theme(legend.justification=c(0,1), legend.position=c(0,1))

Tissue-Wide joint global estimates

tislist <- scan(my.dir %&% 'TW.ten.tissue.list',sep="\n",what="character")
tw <- data.frame()
for(tis in tislist){
  data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")  
  glo.jt <- select(data,tissue,glo.jt.h2,glo.jt.se) %>% arrange(glo.jt.h2) %>% mutate(glo.jt.h2=ifelse(is.na(glo.jt.h2),0,glo.jt.h2)) %>% mutate(ymin = pmax(0, glo.jt.h2 - 2 * glo.jt.se), ymax = pmin(1, glo.jt.h2 + 2 * glo.jt.se)) %>% mutate(nonzeroCI=ymin>0,place=1:length(data$tissue))
  #glo.jt <- glo.jt[complete.cases(glo.jt),]
  tw <- rbind(tw,glo.jt)
}

twrmNA <- tw[complete.cases(tw),]

p<-ggplot(twrmNA,aes(x=place,y=glo.jt.h2,ymin=ymin,ymax=ymax,color=nonzeroCI) ) + facet_wrap(~tissue,ncol=2) + geom_pointrange(col='gray')+geom_point()+ylab(expression("global h"^2)) + xlab(expression("genes ordered by joint global h"^2))+theme_bw() + coord_cartesian(ylim=c(0,1)) + ggtitle("CT and TW")

###calc % nonzero for each tissue
### ADD to plot legend
pvec<-vector()
h2vec<-vector()
a<-twrmNA %>% select(tissue,nonzeroCI) %>% spread(tissue,nonzeroCI)
for(i in 1:10){
  tis<-colnames(a)[i]
  cat("\n\n---",tis,"---\n")
  print(table(a[,i]))
  per <- signif(table(a[,i])/sum(table(a[,i])),3)*100
  print(per)
  pvec <- c(pvec,per[2])
}
## 
## 
## --- cross-tissue ---
## 
## FALSE  TRUE 
## 12123   333 
## 
## FALSE  TRUE 
## 97.30  2.67 
## 
## 
## --- Adipose-Subcutaneous ---
## 
## FALSE  TRUE 
##  3553    16 
## 
##  FALSE   TRUE 
## 99.600  0.448 
## 
## 
## --- Artery-Tibial ---
## 
## FALSE  TRUE 
##  3499     4 
## 
##  FALSE   TRUE 
## 99.900  0.114 
## 
## 
## --- Heart-LeftVentricle ---
## 
## FALSE 
##  2167 
## 
## FALSE 
##   100 
## 
## 
## --- Lung ---
## 
## FALSE  TRUE 
##  4226     3 
## 
##   FALSE    TRUE 
## 99.9000  0.0709 
## 
## 
## --- Muscle-Skeletal ---
## 
## FALSE  TRUE 
##  3912   144 
## 
## FALSE  TRUE 
## 96.40  3.55 
## 
## 
## --- Nerve-Tibial ---
## 
## FALSE  TRUE 
##  3066     1 
## 
##    FALSE     TRUE 
## 100.0000   0.0326 
## 
## 
## --- Skin-SunExposed(Lowerleg) ---
## 
## FALSE  TRUE 
##  3513    26 
## 
##  FALSE   TRUE 
## 99.300  0.735 
## 
## 
## --- Thyroid ---
## 
## FALSE  TRUE 
##  3508     2 
## 
##  FALSE   TRUE 
## 99.900  0.057 
## 
## 
## --- WholeBlood ---
## 
## FALSE  TRUE 
##  3518   153 
## 
## FALSE  TRUE 
## 95.80  4.17
###calc mean h2 for each tissue
a<-twrmNA %>% select(tissue,glo.jt.h2) %>% spread(tissue,glo.jt.h2)
for(i in 1:10){
  tis<-colnames(a)[i]
  meanh2 <- signif(mean(a[,i],na.rm=TRUE),3)
  cat(tis,"mean h2:",meanh2,"\n")
  h2vec <- c(h2vec,meanh2)
}
## cross-tissue mean h2: 0.145 
## Adipose-Subcutaneous mean h2: 0.247 
## Artery-Tibial mean h2: 0.268 
## Heart-LeftVentricle mean h2: 0.343 
## Lung mean h2: 0.318 
## Muscle-Skeletal mean h2: 0.192 
## Nerve-Tibial mean h2: 0.27 
## Skin-SunExposed(Lowerleg) mean h2: 0.223 
## Thyroid mean h2: 0.24 
## WholeBlood mean h2: 0.252
ann_text <- data.frame( glo.jt.h2 = rep(0.78,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2<-p+geom_text(data=ann_text,aes(label=paste("mean_h^2 ==",mean_h2,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3) 
ann_text <- data.frame( glo.jt.h2 = rep(0.9,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2+geom_text(data=ann_text,aes(label=paste("percent_TRUE ==",percent,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)+ theme(legend.justification=c(0,1), legend.position=c(0,1))

Cross-tissue v. Tissue-specific joint h2 and se

tislist <- scan(my.dir %&% 'ten.tissue.list',sep="\n",what="character")
ts <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tislist[1] %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t") 

##LOCAL
tsh2 <-ts %>% select(ensid,loc.jt.h2) 
colnames(tsh2) = c("ensid",tislist[1])
tsse <-ts %>% select(ensid,loc.jt.se) 
colnames(tsse) = c("ensid",tislist[1])

for(tis in tislist[2:10]){
  data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.TS.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")  
  datah2 <-data %>% select(ensid,loc.jt.h2) 
  colnames(datah2) = c("ensid",tis)
  datase <-data %>% select(ensid,loc.jt.se) 
  colnames(datase) = c("ensid",tis)
  tsh2 <- inner_join(tsh2,datah2,by="ensid")
  tsse <- inner_join(tsse,datase,by="ensid")
}

gtsh2<-gather(tsh2,"CrossTissue","Tissue",3:11)
colnames(gtsh2) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtsh2,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red')  + ylab(expression('Cross-Tissue h'^2)) + xlab(expression('Tissue-Specific Joint local h'^2)) + coord_cartesian(ylim=c(-0.1,1.1),xlim=c(-0.1,1.1)) + theme_bw() 

gtsse<-gather(tsse,"CrossTissue","Tissue",3:11)
colnames(gtsse) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtsse,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red')  + ylab('Cross-Tissue SE') + xlab('Tissue-Specific Joint local SE') + coord_cartesian(ylim=c(-0.01,0.18),xlim=c(-0.01,0.18)) + theme_bw() 

##GLOBAL
tsh2 <-ts %>% select(ensid,glo.jt.h2) 
colnames(tsh2) = c("ensid",tislist[1])
tsse <-ts %>% select(ensid,glo.jt.se) 
colnames(tsse) = c("ensid",tislist[1])

for(tis in tislist[2:10]){
  data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.TS.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")  
  datah2 <-data %>% select(ensid,glo.jt.h2) 
  colnames(datah2) = c("ensid",tis)
  datase <-data %>% select(ensid,glo.jt.se) 
  colnames(datase) = c("ensid",tis)
  tsh2 <- inner_join(tsh2,datah2,by="ensid")
  tsse <- inner_join(tsse,datase,by="ensid")
}

gtsh2<-gather(tsh2,"CrossTissue","Tissue",3:11)
colnames(gtsh2) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtsh2,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red')  + ylab(expression('Cross-Tissue h'^2)) + xlab(expression('Tissue-Specific Joint global h'^2)) + coord_cartesian(ylim=c(-0.1,1.1),xlim=c(-0.1,1.1)) + theme_bw() 

gtsse<-gather(tsse,"CrossTissue","Tissue",3:11)
colnames(gtsse) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtsse,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red')  + ylab('Cross-Tissue SE') + xlab('Tissue-Specific Joint global SE') + coord_cartesian(ylim=c(-0.01,1.4),xlim=c(-0.01,1.4)) + theme_bw() 

Cross-tissue v. Tissue-wide joint h2 and se

tislist <- scan(my.dir %&% 'rmTW.ten.tissue.list',sep="\n",what="character")
tw <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tislist[1] %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t") 

##LOCAL
twh2 <-tw %>% select(ensid,loc.jt.h2) 
colnames(twh2) = c("ensid",tislist[1])
twse <-tw %>% select(ensid,loc.jt.se) 
colnames(twse) = c("ensid",tislist[1])

for(tis in tislist[2:10]){
  data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.TW.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")  
  datah2 <-data %>% select(ensid,loc.jt.h2) 
  colnames(datah2) = c("ensid",tis)
  datase <-data %>% select(ensid,loc.jt.se) 
  colnames(datase) = c("ensid",tis)
  twh2 <- inner_join(twh2,datah2,by="ensid")
  twse <- inner_join(twse,datase,by="ensid")
}

gtwh2<-gather(twh2,"CrossTissue","Tissue",3:11)
colnames(gtwh2) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtwh2,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red')  + ylab(expression('Cross-Tissue h'^2)) + xlab(expression('Tissue-Wide Joint local h'^2)) + coord_cartesian(ylim=c(-0.1,1.1),xlim=c(-0.1,1.1)) + theme_bw() 

gtwse<-gather(twse,"CrossTissue","Tissue",3:11)
colnames(gtwse) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtwse,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red')  + ylab('Cross-Tissue SE') + xlab('Tissue-Wide Joint local SE') + coord_cartesian(ylim=c(-0.01,0.21),xlim=c(-0.01,0.21)) + theme_bw() 

##GLOBAL
twh2 <-tw %>% select(ensid,glo.jt.h2) 
colnames(twh2) = c("ensid",tislist[1])
twse <-tw %>% select(ensid,glo.jt.se) 
colnames(twse) = c("ensid",tislist[1])

for(tis in tislist[2:10]){
  data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.TW.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")  
  datah2 <-data %>% select(ensid,glo.jt.h2) 
  colnames(datah2) = c("ensid",tis)
  datase <-data %>% select(ensid,glo.jt.se) 
  colnames(datase) = c("ensid",tis)
  twh2 <- inner_join(twh2,datah2,by="ensid")
  twse <- inner_join(twse,datase,by="ensid")
}

gtwh2<-gather(twh2,"CrossTissue","Tissue",3:11)
colnames(gtwh2) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtwh2,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red')  + ylab(expression('Cross-Tissue h'^2)) + xlab(expression('Tissue-Wide Joint global h'^2)) + coord_cartesian(ylim=c(-0.1,1.1),xlim=c(-0.1,1.1)) + theme_bw() 

gtwse<-gather(twse,"CrossTissue","Tissue",3:11)
colnames(gtwse) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtwse,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red')  + ylab('Cross-Tissue SE') + xlab('Tissue-Wide Joint global SE') + coord_cartesian(ylim=c(-0.01,1.6),xlim=c(-0.01,1.6)) + theme_bw()